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Abstract 

We provide a proof that when there is no forcing the recently introduced projective 
dynamics Monte Carlo algorithm gives the exact lifetime of the metastable state, 
within statistical uncertainties. We also show numerical evidence illustrating that 
for slow forcing the approach to the zero-forcing limit is rather rapid. The model 
studied numerically is the 3-dimensional 3-state Potts ferromagnet. 
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One of the most challenging hurdles to be overcome in computational mate- 
rials science is finding a direct way to bridge the disparate time and length 
scales required in simulations of models of materials. The length scales go 
from the atomic to the size of the actual material sample. The shortest time 
scale in many instances can be associated with an inverse phonon frequency, 
and the longest time scale with the mean time before macroscopic material 
failure. One way of bridging these time and length scales is to devise better al- 
gorithms. Such algorithms need to increase computational efficiency by orders 
of magnitude in order to allow the desired material simulations to proceed on 
today's computers. 

Although there are a number of efficient algorithms for the statics of model 
systems [1], the number of advanced algorithms for dynamics is more limited. 
These include the ra-fold way algorithm [2], and its generalization, the Monte 
Carlo with Absorbing Markov Chains (MCAMC) algorithm [3] . 

Recently, the authors have devised an improved dynamical algorithm for mod- 
els with discrete numbers of states at each lattice site. This algorithm is called 
the Projective Dynamics (PD) method [4]. The PD algorithm builds on previ- 
ous work using projections onto slow variables in dynamical models [5,6]. The 
PD algorithm has been shown to allow simulations of the dynamics of systems 
such as the Ising model for very long times and large lattice sizes [7] , requiring 
quite modest computer resources. It has also been shown to allow identification 
of physical quantities, such as activation volumes and Barkhausen volumes, in 
dynamical models of magnetic domain growth [8]. 

In this brief paper we concentrate on the use of the PD algorithm to obtain 
the average lifetime of the metastable state in a model system. In particular, 
we apply the PD algorithm to the escape from the metastable state of the 
three-dimensional 3-state Potts ferromagnet. 

The 3-state Potts model in applied field H has the Hamiltonian 

n = -J E °i) + H E W, °i) - (!) 

where <7j G {0,1,2} is the 'spin' at lattice site %. The first summation runs 
over all nearest-neighbor pairs on a simple-cubic lattice, and J>0 corresponds 
to ferromagnetic interactions. The second sum runs over all lattice sites, and 
it only changes the energies of spins in the states and 1. The system can be 
characterized by the concentrations {n ,ni,n 2 }, J2 n i=l, of spins in the three 
states. 

The Monte Carlo simulation uses local updates at randomly chosen sites. We 
use simple cubic lattices of size L 3 =V with periodic boundary conditions. 
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In this paper we use Glauber dynamics, and we measure the time in Monte 
Carlo Steps per Spin (MCSS). The simulation starts with all spins in state 0, 
with the temperature T below the critical temperature T c . We measure the 
average lifetime r, defined as the average first-passage time to a configuration 
with half of the spins in the stable state, which here is state 1. Our Monte 
Carlo simulation can be viewed as a random walk which starts in the state 
with riQ=V and goes until it reaches a state with n—N—V/2, where we have 
dropped the subscript on n=n 1 . 

The spins in the current configuration can be divided into classes specified by 
the state a of the spin and by the numbers {a, b, 6 — a — b} of neighbors in 
the states {0,1,2}. Let c^ b (rio, rii, n 2 ) be the average concentration of spins 
in the class {a, a, b}, conditional on the total concentrations rii of spins in 
state i. This average is obtained from repeated metastable simulations. The 
class concentrations are normalized such that J2aab c ab( n 0i n 2)=l for all 
{n , rii, n 2 }. Let p° b be the probability that a spin in the class {a, a, b} will flip 
to the state a' when visited by the updating algorithm. We project the spin 
class concentrations onto bin n by defining c° b (n)=Y /no ,n 2 c ab( n 0i n i=n, n 2 ). 
The projected flip rates are defined in terms of the c° b (n) and the flipping 
probabilities p°f as: 

9{n) = £ < fc (n) pll , s{n) = £ c\ b {n) p 1 *. (2) 

ab,ijytl ab,a^l 



The rates s(n) and g(n) correspond to the shrinkage and growth rates of the 
stable phase, respectively. They depend on how the configurations are gener- 
ated to obtain the c^ b {n). Here we obtain the configurations during simulation 
of escape from the metastable state, with and without the forcing described 
below. 

Using standard methods from the theory of absorbing Markov chains [3,6], 
the mean lifetime r and the total average time h(n) (measured in MCSS, with 
h(N) = 0) spent by the random walker in the state n are given, as shown in 
the appendix, by 

r =|>(«), Hn-l)= g(n _ 1) . (3) 



An enhancement to the PD algorithm described above was presented in Ref . [4] . 
At any time, the system is only allowed to have the number of spins in the sta- 
ble phase, n, strictly larger than a lower bound n min (t). The simulation starts 
with n min =— 1, and n m j n is increased at a slow constant rate. The class concen- 
trations c° b {n) are measured during this forced escape and used to calculate 
r from Eq. (3). In the appendix we prove that if the random walk does not 
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hit the moving lower bound, the value of r obtained with the PD algorithm is 
exact. Figure 1 shows the approach of r to the slow- forcing limit as the forc- 
ing rates is lowered. Even at moderately high forcing rates the lifetime (with 
statistics between 10 2 and 10 3 escapes from the metastable state) is indistin- 
guishable from the value at zero forcing rate. The forcing, however, decreases 
substantially the amount of computer time required to obtain r. In the con- 
ventional direct simulation approach, the necessary computer time increases 
exponentially with the depth of the metastable free-energy minimum. On the 
other hand, in the forced-escape method, the computer time is proportional 
to the depth of the free-energy minimum itself. In other words, the proposed 
method allows calculations of an exponentially hard problem in linear time. 
Although here we only present data for a temperature and field at which di- 
rect simulation is possible and the speed-up is only about 10 2 -10 3 , the forced 
escape method allows one to obtain lifetimes many orders of magnitude longer 
than would be possible without forcing [4,7]. 

In conclusion, we have obtained the average lifetime, r, of the metastable 
state of the 3-state three-dimensional Potts ferromagnet using forcing in the 
PD (Projective Dynamics) method. In the limit of zero forcing, we prove 
in the appendix that the value of r obtained is exact. Substantial savings in 
computer time can be achieved using relatively moderate forcing rates without 
significantly affecting the obtained values of r. 
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APPENDIX: Lifetimes from projective dynamics 

We prove that if there is no forcing, the lifetime obtained from the projective 
dynamics algorithm is exact. Of course, there will still be statistical uncertain- 
ties in the lifetime obtained from any finite length simulation. The proof is 
given for the model under consideration in this paper, the three-dimensional 
three-state Potts model. Furthermore, we choose state 1 as the slow variable 
on which we project the dynamics. The proof can easily be generalized to 
other models and dimensions as long as the number of states in the model is 
finite. The notation we use is given below. 

• V — total number of spins in the system. 

• S — \Jn=o S n — partition of the state space S of the model into disjoint 
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subspaces where S n contains all microstates with exactly n spins in state 1. 

• [n, i] - - label for the ith state in S n , i = 1, 2, . . . , \S n \, where | • • • | denotes 
the number of elements in the partition. 

• Cab( n i'i) — the number of those spins of the microstate [n, i] which are in 
state a and have a and b neighbors in states and 1, respectively. Thus, 
EabaCabi^i) = V for each [n,i]. 

• (% b (n, i) = C° b (n,i)/V denotes concentrations. 

• Vab ~ probability that a chosen spin in state a having a, b neighbors in 
states and 1 will flip to state a'. 

Consider direct Monte Carlo measurement of the mean lifetime r without 
forcing. This can be accomplished by setting n m i n = — 1 at all times. Start 
an escape from the state [0, 1] and continue the simulation until the Monte 
Carlo "random walker" reaches one of the absorbing states [n st0 p,j] with 
j—1,2, . . . , (S'nstopl- A total number iV esc of independent escapes are realized. 
Denote by W(n,i) the number of visits in the state [n,i] generated by this 
repeated measurement (counting rejected Monte Carlo moves as new visits in 
the same state). Then 

1 \Sn\ 

h(n) = -J2W(n,t)/N esc (4) 



is the average (sojourn) time spent in S n , in units of Monte Carlo Steps per 
Spin (MCSS). The mean lifetime is obtained by 

^stop 1 -i ^stop 1 |£Vi| 

r= E Kn) = V/ E EWM/N esc . (5) 

n=0 V n=0 i=l 

One can imagine that each random walker generates an oriented "world line" 
starting in [0, 1] and ending in some [n stop , j], because the transitions are possi- 
ble only between neighboring S n and subspaces. Therefore, an "equation 
of continuity" holds for each n < n stop : 

-^n+l — N csc + N n+ i^ n (6) 

where stands for the number of transitions between the subspaces with 

t and £' spins in state 1. Equation (6) can be written in terms of W, C° b and 
PaZ . To do this we observe that the probability that a particular attempted 
spin flip results in a move from S n to S n+ i is V^ 1 J2 a b C° b (n, i)p°l, where 
the factor V^ 1 represents the probability of choosing one particular spin. Con- 
sequently, 
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^f>(MEE^(M^ = 

V i=l ab ct^I 

1 l^n+ll 

a^csc + 7 F E H '( n +M)EE + !> Op£ (7) 

K i=l a& cr^l 

Since each must be visited at least iV esc times, we can multiply the two 
terms involving sums by an identity and write the above equation as 



i=i 



In Eq. (8) one recognizes the average lumped class concentration 



< b (n) = E W(n, iK b (n, i) /]T W(n, k) (9) 
i=i / fe=i 



which enters the definition of the global spin-flip rates, Eq. (2). We divide by 
VN esc to express the sojourn times in units of MCSS. This gives 

Kn) E E PlKbin) = V- 1 +% + l)EE Plbclbin + 1), (10) 

ab <r^l ab cr^l 



which leads directly to Eq. (3). Thus, in the absence of forcing, our formula 
for the lifetime becomes exact. 
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Figure Captions 



Fig. 1. The measured average lifetime, r, is shown as a function of the inverse 
forcing rate. This figure illustrates that excellent values for r are obtained 
even with moderately fast forcing. The figure is for 10 2 -10 3 escapes from the 
metastable state. To achieve uniform accuracy, more escapes are needed at 
high forcing rates. 
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